The Fragmentation of Cores and the Initial Binary Population 

Simon P Goodwin 

University of Sheffield 

Pavel Kroupa 

Rheinische Fredrich-Wilhelms-Universitat Bonn 

Alyssa Goodman 

Harvard-Smithsonian Center for Astrophysics 

Andreas Burkert 

Ludwig-Maximilhan-Universitat Munich 

Almost all young stars are found in multiple systems. This suggests that protostellar cores 
almost always fragment into multiple objects. The observed properties of multiple systems 
such as their separation distribution and mass ratios provide strong constraints on star formation 
theories. We review the observed properties of young and old multiple systems and find that 
the multiplicity of stars changes. Such an evolution is probably due to (a) the dynamical decay 
of small-A'^ systems and/or (b) the destruction of multiple systems within dense clusters. We 
review simulations of the fragmentation of rotating and turbulent molecular cores. Such models 
almost always produce multiple systems, however the properties of those systems do not match 
observations at all well. Magnetic fields appear to supress fragmentation, prehaps suggesting 
that they are not dynamically important in the formation of multiple systems. We finish by 
discussing possible reasons why theory fails to match observation, and the future prospects for 
this field. 



1. INTRODUCTION. 

CoiTectly predicting the properties of young multiple 
systems is one of the most challenging tests of any theory 
of star formation. In this chapter we discuss the current 
understanding of how dense prestellar cores fragment into 
multiple stars including brown dwarfs (by 'stars' we gener- 
ally mean both stars and brown dwarfs). 

Firstly, we will discuss the important observed proper- 
ties of both young and old binary systems and the differ- 
ences between them. Then we will describe the possible 
origins of the differences between the young and old sys- 
tems and hopefully convince the reader that almost all stars 
must form in multiple systems. This initial binary popu- 
lation must form from the fragmentation of star-forming 
dense molecular cores and so we discuss the observed prop- 
erties of these cores that may influence their ability to form 
multiple systems. We will then review the current mod- 
els of core fragmentation, with an emphasis on turbulence 
as the mechanism that promotes fragmentation. Finally we 
will examine why theory currently fails to coiTectly predict 
binary properties. 

2. THE PROPERTIES OF MULTIPLE SYSTEMS. 

There has been an extensive study of binary properties 
over the past two decades with the modern study often 



marked as beginning with the detailed survey by Duquen- 
noy and Mayor (1991, hereafter DM91). Multiple systems 
in the field are by far the best studied due to the availabil- 
ity of local samples whose completeness is easier to esti- 
mate, and the properties of the field provide the benchmark 
against which younger samples are measured. 

2.1 Multiple Systems in the field. 

2.1.1 Multiplicity fraction. The fraction of field stars in 
multiple systems is found to be high and increases with in- 
creasing primary mass. Many different measures are used 
to quantify multiplicity which can become very confus- 
ing. Most important is the 'multiplicity frequency' /mult = 
(Bh-Th-Qh-...)/(Sh-Bh-Th-Qh-...) where S, B, T and Q are the 
numbers of single, binary, triple and quadruple systems re- 
spectively, thus /mult represents the probability that any 
system is a multiple system (see Reipurth and Zinnecker, 
1993). 

The raw value of the field G-dwarf multiplicity fre- 
quency found by DM91 is 0.49, when corrected for incom- 
pleteness this rises to 0.58. However, recent studies using 
Hipparcos data have shown that it may even be higher than 
this (Quist and Lindegren, 2000; Sdderhjelm, 2000). 

It should be noted that the brown dwarf (BD) binary frac- 
tion has generally been considered to be much lower than 
that of stars at ~ 0.1 — 0.2 (Bouy et ai, 2003; Close et at.. 
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2003; Gizis et al, 2003; Martin et al, 2003). However, 
recent studies have suggested that the BD muhipHcity fre- 
quency may be significantly higher, possibly exceeding 0.5 
(Pinfield et al., 2003; Maxted and Jeffries, 2005), provided 
most BDs reside in very tight BD-BD pairs. 

2.7.2 Separation distribution. The binary separation dis- 
tribution is very wide and flat, usually modelled as a log- 
normal with mean ^ 30 AU and variance criogd ^ 1-5 
(DM91 for G-dwarfs). This is illustrated in Fig. 1 where 
the field period distribution is compared to that of young 



stars (note that cjr/P 



where a is in AU, P is in 



years and mgys is the system mass in Mq). 

A similar distribution is found for M-dwarfs by Fischer 
andMarcy (1992), and generally seems to hold for all stars, 
although the maximum separations do appear to decrease 
somewhat, but not substantially, for stars with decreasing 
mass (Close et al, 2003). Very low-mass stars (VLMSs) 
and BDs seem strongly biased towards very close compan- 
ions with semi-major axes a < a,„ax ~ 15 AU (Close et 
al., 2003; Gizis et at., 2003; Pinfield et al., 2003; Maxted 
and Jeffries, 2005) in contrast to those of stars that have 
flmax <; 100 AU (DM91; Fischer and Marcy, 1992; Mayor 
et al, 1992). It is this unusual separation distribution which 
may have led to the underestimate of the BD multiplicity 
fraction. The much smaller Omax for VLMSs and BDs com- 
pared to the other stars cannot be a result of disruption in a 
cluster environment but must be due to the inherent physics 
of their formation (Kroupa et al., 2003). 

2.1.3 Mass ratio distribution. DM9 1 found that Galactic- 
field systems with a G-dwarf primary have a mass-ratio 
distribution biased towards small values such that it does 
not follow the stellar IMF which would predict a far 
larger number of companions with masses m2^0.3MQ 
(Kroupa, 1995b). For short-period binaries, the mass-ratio 
distribution is biased towards similar-mass pairs (Mazeh 
et ah, 1992). Integrating over all periods, for a sam- 
ple of nearby systems with primary masses in the range 
0.1;$mi/MQ <1, Reid and Gizis (1997) find the mass- 
ratio distribution to be approximately flat and consistent 
with the IMF (Fig. 2). 

2.7.4 Eccentricity distribution. Binary systems have a 
thermalised eccentricity distribution (Eqn. 3 below) for pe- 
riods P ^ 10^~* d, with tidally circularised binaries dom- 
inating at low separations (DM91; Fischer and Marcy, 
1992). 

2.1.5 Higher-order systems. DM91 find the uncorrected 
ratio of systems of different multiphcity to be S:B:T:Q = 
1.28 : 1 : 0.175 : 0.05 (see also Tokovinin and Smekhov, 
2002), suggesting that roughly 20% of multiple systems are 
high-order systems. Concerning the origin of high-order 
multiple systems in the Galactic field, we note that many, 
and perhaps most, of these may be the renmants of star clus- 
ters (Goodwin and Kroupa, 2005). 

2.2 Pre-Main Sequence Multiple Systems. 

The properties of pre-main sequence (PMS) multiple 
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Fig. 1. — The period distribution function. Letters show the 
observed fraction of field G, K and M-stars and PMS stars (P). 
The solid curve shows the model initial period distribution (see 
Eqn. 2). The light histogram is the initial binary population in the 
simulations of Kroupa (1995b) which evolves through dynamical 
interactions in a cluster into a field-like distribution shown by the 
heavy histogram. 



systems are much harder to determine than those in the 
field. We refer the reader to the chapter by Duchene et al. 
for a detailed review of the observations of PMS multiple 
systems and the inherent problems. 

Probably the most important difference between the 
PMS and field populations is that young stars have a signif- 
icantly higher multiplicity fraction than the field (see the 
chapter by Duchene et al.; also see Fig. 1). 

The separation distribution of PMS stars also appears 
different to that in the field with an over-abundance of bina- 
ries with separations of a few hundred AU (Mathieu, 1994; 
Patience et al., 2002; Fig. 1). More specifically, the binary 
frequency in the separation range 100 — 1000 AU is a fac- 
tor of ~ 2 higher than in the field (Mathieu, 1994; Patience 
et al, 2002; Duchene et al, 2004). Extrapolating this in- 
crease across the whole separation range implies that /mult 
for PMS stars could be as high as 100%. (It appears that in 
Taurus the binary frequency is ~ 100% for stars > O.SM©, 
Leinert et al., 1993; Kdhler and Leinert, 1998). 

The mass-ratio distribution of PMS stars is similar to 
the field population. A detailed comparison is not yet pos- 
sible because low-mass companions to pre-main sequence 
primaries are very difficult to observe as the available re- 
sults depend mostly on direct imaging or speckle interfer- 
ometry, while for main-sequence systems radial-velocity 
surveys have been done over decades (DM91). Thus, using 
near-infrared speckle interferometry observations to obtain 
resolved JHK-photometry for the components of 58 young 
binary systems, Woitas et al. (2001) found that the mass- 
ratio distribution is flat for mass ratios q > 0.2 which is 
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consistent with random pairing from the IMF, i.e. fragmen- 
tation processes rather than common-accretion (Fig. 2). 

The eccentricity distribution of PMS stars also is simi- 
lar to the field, with a thermaUsed distribution except at low 
separations where tidal circularisation has occurred rapidly 
(e.g., Kroupa, 1995b; White and Ghez, 2001). Finally, the 
ages of components in young multiple systems appear to be 
very similar (White and Ghez, 2001). 

It is currently unclear what the proportion of higher- 
order multiples is in young systems (see the chapter by 
Duchene et al). We will revisit this question in the next 
section. 

2.3 The evolution of binary properties. 

The observations described in Sections 2.1 and 2.2 
clearly show that at least the binary fraction and separa- 
tion distributions evolve significantly between young stellar 
populations and the field. 

Indeed, binary properties are seen to change even within 
populations in star forming regions. The binary frac- 
tion is found to vary between embedded and (older) non- 
embedded sources in both Taurus and p Oph (Duchene et 
al, 2004; Haisch et al, 2004). Also, the mass ratio dis- 
tributions of massive stars appear to depend on the age of 
the cluster with those in young clusters being consistent 
with random sampling from the IMF and those in dynam- 
ically evolved populations favouring equal-mass compan- 
ions (Section 2.4). 

The evolution of binary properties has been ascribed 
to two mechanisms. Firstly, the rapid dynamical de- 
cay of young small- A'' clusters within cluster cores (e.g., 
Reipurth and Clarke, 2001; Sterzik and Durisen, 1998, 
2003; Durisen et al, 2001; Hubber and Whitworth, 2005; 
Goodwin and Kroupa, 2005; Umbreit et al, 2005), and 
secondly, the dynamical destruction of multiples by inter- 
actions in a clustered environment can modify an initial 
PMS-like distribution into a field-like distribution (Kroupa, 
1995a, b; Kroupa et al, 2003; Figs. 1, 3). 

2.3.1 Small- N decay. Multiple systems containing N > 
3 stars are unstable to dynamical decay unless they form 
in a strongly hierarchical configuration (stability criteria for 
N > 2 systems are provided by Eggleton and Kiseleva, 
1995). Generally, a triple system is unstable to decay with 
a half-life of 

where R is the size of the system, and Mgtars is the mass 
of the components (Anosova, 1986). The decay time for 
R = 250 AU and Mstars = IM© is - 55 kyr which is of 
order the duration of the embedded phases of young stars, 
thus ejections should mainly occur during the main Class 
accretion phase of PMS objects (e.g., Reipurth, 2000). 
Indeed, one such early dynamical decay appears to have 
been observed by Gomez et al. (2006), and this is probably 



the process at work to reduce the binary fraction between 

the embedded and non-embedded stars seen by Duchene et 
al (2004) and Haisch et al. (2004). These early dynami- 
cal processes cause embedded protostars to be ejected from 
their natal envelopes, possibly causing abrupt transitions of 
objects from class O/I to class II/III (Reipurth, 2000; Good- 
win and Kroupa, 2005). 

Significant numbers of small-iV decays will dilute any 
initial high multiplicity fraction very rapidly to a small bi- 
nary fraction for the whole population (e.g., N = 5 systems 
would lead to a population with a binary fraction / = 1/5 
within < 10^ yr). The observed high binary fraction in 
about 1 Myr old populations thus suggests that the forma- 
tion of TV > 2 systems is the exception rather than the 
rule (Goodwin and Kroupa, 2005). Ejections would occur 
mostly during the very early Class stage such that ejected 
embryos later appear as free-floating single very low-mass 
stars and BDs (Reipurth and Clarke, 2001). However, the 
small ratio of the number of BDs per star, w 0.25 (Munech 
et al, 2002; Kroupa et al, 2003; Kroupa and Bouvier, 
2003b; Luhman, 2004), again suggests this not to be a very 
common process even if all BDs form from ejections. 

Ejections have two main consequences: a significant 
reduction in the semi-major axis of the remaining stars 
(Anosova, 1986; Reipurth, 2000; Umbreit et al, 2005); 
and the preferential ejection of the lowest mass component 
(Anosova, 1986; Sterzik and Durisen, 2003). The early 
ejection of the lowest-mass component forms the basis of 
the embryo ejection scenario of BD formation (Reipurth 
and Clarke, 2001; Bate et al, 2002). 

The A'^-body statistics of the decay of small- systems 
has been studied by a number of authors (Anosova, 1986; 
Sterzik and Durisen, 1998, 2003; Durisen et al, 2001; 
Goodwin et al, 2005; Hubber and Whitworth, 2005). How- 
ever, only Umbreit et al (2005) have attempted to include 
the effects of accretion on the A'^-body dynamics which ap- 
pear to have a significant effect - especially on the degree 
of hardening of the binary after ejection. Goodwin et al 
(2004a, b) and Delgado Donate et al (2004a, b) have sim- 
ulated ensembles of cores including the full hydrodynam- 
ics of star formation, however proper statistical conclusions 
about the effects of ejections are difficult to draw due to the 
different numbers of stars forming in each ensemble (which 
there is no way of controlling a priori), and the smaller num- 
ber of ensembles that may be run in a fully hydrodynamic 
context. However, some conclusions appear from these and 
other studies (Whitworth et al, 1995; Bate and Bonnell, 
1997; Bate et al, 2003; Delgado Donate et al, 2003). 
Firstly, that early ejections are very effective at hardening 
the remaining stars (c.f. Umbreit et al, 2005). Secondly, 
this early hardening tends to push the mass ratios of close 
binaries towards unity. This occurs as the low-mass com- 
ponent has a higher specific angular momentum than the 
primary and so is more able to accrete mass from the high 
angular momentum circumstellar material (see Whitworth 
et al, 1995; Bate and Bonnell, 1997). However, Ochi et al 
(2005) find in detailed 2D simulations of accretion onto bi- 
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naries that the gas accretes mainly onto the primary due to 

shocks removing angular momentum. 

One significant caveat to the previous discussion is that 
the gradual formation (over ~ 0.1 Myr) of stars allows 
far more stable triples and higher-order multiples to form 
than expected observationally, which are stable for at least 
10 Myr (Delgado Donate et al, 2004a, b). Indeed, simu- 
lations that form a large number of stars often form very 
hierarchical higher-order multiples (often quadruples and 
quintuples formed when even larger systems decay) which 
are not observed (Delgado Donate et al., 2003, 2004a, b; 
Goodwin et al, 2004a, b). Such systems would probably be 
destroyed during the cluster destruction phase (see below), 
but not dilute the binary fraction on very short timescales. 

2.3.2 Dynamical destruction in clusters. In the highly- 
clustered environments in which most stars are thought to 
form (e.g., Lada and Lada, 2003) dynamical interactions 
will be common and may disrupt many initially binary sys- 
tems. 

Binaries can be sub-divided into three dynamical groups: 
(i) the wide, or soft, binaries, (ii) the dynamically active 
binaries, and (iii) the tight or hard binaries. 

Wide binaries have orbital velocities much smaller than 
the velocity dispersion, a, in a cluster and are easily dis- 
rupted. This is best seen by a gedanken experiment, where 
we construct a reduced-mass particle (a test particle orbiting 
in a fixed potential with total mass, eccentricity and orbital 
period equal to that of the binary in question) in a heat bath 
of field particles (the cluster stars). If the orbital velocity of 
the test particle is smaller than the typical velocities of the 
field particles {vorh ^ f). then the test particle will gain 
kinetic energy by encounters, i.e. by the principle of energy 
equipartition, until it's orbital velocity surpasses the bind- 
ing energy of the binary. The binary consequently gets dis- 
rupted. Energy conservation requires the heat bath to cool; 
cluster cooUng has been seen in N-hoAy computations by 
Kroupa et al., (1999), but the effect is not significant for 
cluster dynamics. The general effect of this process is that 
binaries with weak binding energies are disrupted (i.e., bi- 
naries with long periods and/or small mass-ratios). 

Hard binaries, on the other hand, can be represented 
by a reduced-mass binary in which the test particle has 
Worb ^ o", so that energy equipartition leads to a reduc- 
tion of Uorb and to an increase of a (cluster heating). This 
increases the binding energy of the binary which heats up 
further (uorb increases as the test particle falls towards the 
potential minimum). This run-away process only stops be- 
cause either the binary merges when it is so tight that the 
constituent stars touch (forming a blue straggler), or be- 
cause the hardened binary receives a re-coU expelUng it 
from the cluster, or the hardening binary evolves to a cross- 
section so small that the binary becomes essentially unre- 
solved in further interactions. 

Heggie (1975) and Hills (1975) studied the details of 
these processes and formulated the Heggie-Hills law of bi- 
nary evolution in clusters: "soft binaries soften and hard 
binaries harden". An important implication of this law is 



that hard binaries can absorb the entire binding energy of 
a cluster and drive the evolution of the core of a massive 
cluster. 

Not accessible to analytical work are active binaries 
with intermediate binding energies. Only full-scale iV-body 
computations can treat the dynamics of the interactions ac- 
curately (e.g., Heggie et al., 1996). Such binaries couple 
efficiently to the cluster, and efficiently exchange energy 
with it. The binary-binary and binary-single-star interac- 
tions form complex resonances and short-Uved higher-order 
configurations that decay by expelling typically the least 
massive member. Active binaries are thus quite efficient 
in exchanging partners, but more work needs to be done 
in order to quantify the exchange rates for typical Galactic 
star-forming clusters. 

The dynamical interactions between binaries or single 
stars will continue to alter the binary properties of the pop- 
ulation as long as it remains relatively dense (i.e. until the 
cluster dissolves or expands significantly after residual gas 
removal). In a series of papers (Kroupa, 1995a, b; Kroupa 
et al., 2001) it has been shown that a population initially 
composed entirely of binary systems with a PMS binary 
separation distribution can evolve into the field-like distri- 
bution through dynamical encounters, as is exemplified by 
the evolution of the mass-ratio distribution (Fig. 2). 

Importantly, however, dynamical interactions in a clus- 
ter cannot form a significant number of binaries from an 
initially single star population, or widen an initially narrow 
separation distribution (Kroupa and Burke rt, 2001). Also, 
the clusters within which most Galactic-field stars form do 
not sufficiently harden an initially wide separation distribu- 
tion to be consistent with the number of tight binaries. Such 
clusters that would lead to significant hardening would dis- 
rupt too many of the wide binaries diluting the Galactic- 
field binary-star population. 

This indicates that the observed broad period distribu- 
tion (Fig. 1) is already imprinted at the time of binary for- 
mation. The existing Af-body simulations of young clusters 
tend to assume a relatively well-mixed and relaxed initial 
distribution of stars. However, real young clusters tend to 
be lumpy and unrelaxed which alters the binary-binary in- 
teraction rate (e.g., Goodwin and Whitworth, 2004), but is 
not expected to change the general outcome (a smaller- AT 
system with a smaller radius is dynamically equivalent to a 
larger- Af system with a larger radius, (Kroupa, 1995a). 

2.3.3 The effect of dynamics on the binary population. 
Small- Af decay will act to modify the birth binary popula- 
tion (ie. that produced by star formation) by reducing the 
overall binary fraction and by hardening the remaining bi- 
naries on a timescale of < 10^ yrs. Small- A^ decay cannot 
occur very often as it would produce too many single PMS 
stars and too many hard binaries (and quite possibly too 
many BDs). However, it is unclear if small- A^ decay is rare 
because cores usually form only binaries and only some- 
times triples, or because higher order multiples are formed 
in stable, hierarchical systems. 

The initial binary population (ie. the population af- 
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Fig. 2. — The mass-ratio distribution of all late-type primary 
stars. The solid dots are observational data by Reid and Gizis 
(1997), while the expected initial distribution is shown as the 
dashed histogram. In a typical star cluster it evolves to the solid 
histogram which reproduces the observed data quite well (from 
Kroupa et al, 2003). 

ter small- A/^ decay and internal energy re-distribution pro- 
cesses, "eigenevolution"-see below, has modified it) is fur- 
ther changed by dynamical interactions within a cluster on a 
timescale of a Myr for typical Galactic star-forming clusters 
{Kroupa, 1995a, b; Adams and Myers, 2001). Encounters 
will destroy soft and active binaries leaving mostly the hard 
binary population unchanged. Crucially it cannot produce 
more binaries in any significant numbers. 

Thus both of the processes that act to modify the birth 
binary population into the field binary population reduce 
the binary fraction. We are led to the conclusion that the 
birth binary fraction must be higher than that of the field. 
This conclusions agrees well with observations. 

Goodwin and Kroupa (2005) argue that the birth binary 
fraction must be almost unity for all stars, and that the low 
binary fraction amongst M-dwarfs is due to the preferen- 
tial destruction of low-mass binaries. Lada (2006), how- 
ever, argues that most stars form as single stars, as most 
M dwarfs are single and most stars are M dwarfs. The cru- 
cial issue here is how many initial M dwarf binaries decay? 
We note that in a model that assumes that stars are born 
with a 60 % binary fraction with companion masses se- 
lected randomly from the IMF and without dynamical dis- 
solution of the binaries leads to a population with the ob- 
served binary fraction- spectral type relationship {Kroupa et 
al, 1993, their Fig. 1 1). In models that assume a 100 % ini- 
tial binary fraction, processing through an 'average' cluster 
{Lada and Lada, 2003) also converts the binary population 
into the observed field population. Observations of the bi- 
narity of a very young (embedded?) M star population in an 
average cluster (as opposed to Taurus which is atypical) are 



Fig. 3. — Evolution of the total binary fraction as a function of 
time for four cluster models with decreasing half-mass radii 7?o.5 
from top to bottom. Each cluster initially contains 200 late-t5'pe 
binaries (with periods, eccentricities and mass-ratios consistent 
with the initial distribution functions described in Section 2.4) ex- 
cept the bottom simulation which contains no binaries. The area 
marked 'GF' shows the observed field binary fraction (DM91). 
The cluster with iio.s = 0.8 pc evolves into a field-like distri- 
bution. With no initial binaries capture is unable to produce a 
significant number of binaries. From Kroupa (1995a). 

required to fully resolve this issue. However, for higher- 
mass stars, it seems almost certain that the initial binary 
fraction must be significantly higher than the field binary 
fraction due to their current relatively high binary fraction. 
It is worth noting that the observed population of high-order 
multiples may be formed from the final handful of stars that 
remain at the end of cluster dissolution which tend to be in 
high-order heiarchical systems (see Goodwin and Kroupa, 
2005). 

Thus, whilst more work needs to be done, especially on 
the effects of the cluster density and on the initial distri- 
bution in star-forming regions with very different physical 
properties, a picture has emerged in which the observed 
high multiplicity PMS population can be modified by sec- 
ular dynainical evolution to produce the field population. 
Most (especially K-dwarfs and later) stars therefore form in 
multiple - probably binary and triple - systems with a very 
wide separation distribution and a relatively flat mass-ratio 
distribution. 

2.4 The initial binary population. 

Given these results, a useful working hypothesis appears 
to be that the initial binary- star properties are invariant to 
star-formation conditions. The observed differences be- 
tween binary populations result from different secular dy- 
namical histories of the respective populations: ie. due to 
different cluster masses and densities. 
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Kroupa (1995b) suggests that the field binary properties 
can be understood if the birth binary population has the fol- 
lowing semi-empirical distribution functions: 
(1 .) Companion masses are chosen randomly from the IMF; 
(2.) The distribution of periods is independent of primary 
mass, and can be described with the following functional 
form, 

IP -1 

where dNip = A^tot fip dlP is the number of binaries with 
periods in the range IP, IP + dlP (IP = log^o-P. Pin days) 
and A^tot is the number of single-star and binary-star sys- 
tems in the sample under consideration; 
(3.) The eccentricity distribution is thermal (all binding en- 
ergies are equally occupied), 

dN = fof{e)de = fo2ede, (3) 

being the number of binary systems with eccentricity in the 
range e, e + de. 

These birth distributions need to be modified for short- 
period binaries (IP ^ 3) through the evolution of the bind- 
ing energy and angular momentum owing to dissipative 
processes within the young binary system termed collec- 
tively as pre-main sequence eigenevolution. This then gives 
the initial distributions which are evident in dynamically 
unevolved populations (e.g., Taurus, Kroupa and Bouvier, 
2003a) and can be used as the initial binary-star population 
in A^-body modelhng of stellar populations. 

The distribution over binding energies and specific- 
angular momenta can be evaluated readily given the above 
distribution functions. Fig. 4 shows that the distribution of 
specific angular momenta of molecular cloud cores forms 
a natural extension to the distribution of specific angular 
momenta of the initial binary stellar population, possibly 
suggesting an evolutionary connection. 

These semi-empirical distribution functions have been 
formulated for late-type stars (primary mass rrip ^ 1 Mq) 
as it is for these that we have the best observational con- 
straints. It is not clear yet if they are also applicable to 
massive binaries. Baines et al. (2006) report a very high 
(/ « 0.7 ± 0.1) binary fraction among Herbig Ae/Be stars 
with the binary fraction increasing with increasing primary 
mass. Furthermore, they find that the circumbinary discs 
and the companions appear to be co-planar thereby support- 
ing a fragmentation origin rather than collisions or capture 
as the origin of massive binaries. Most O stars are believed 
to exist as short-period binaries with q k. 1 (Garcia and 
Mermilliod, 2001), at least in rich clusters, while small q 
appear to be favoured in less substantial clusters such as the 
Orion Nebula Cluster (ONC), being consistent there with 
random pairing (Preibisch et ai, 1999). Kouwenhoven et 
al. (2005) report the A and late-type B binaries in the Scor- 
pius OB2 association to have a mass-ratio distribution not 
consistent with random pairing. The lower limit on the bi- 
nary fraction is 0.52. Perhaps the massive binaries in the 
ONC represent the primordial population, whereas in rich 
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Fig. 4. — The observed specific angular momentum distribution 
of molecular cloud cores by Goodman et al. (1993) is shown as 
the rightmost dashed histogram, while the initial binary-star pop- 
ulation (Section 2.4) is plotted as the solid histogram. It evolves 
to the dot-dashed histogram after passing through a typical star 
cluster (from Kroupa, 1995b). 

clusters and in OB associations the population has already 
dynamically evolved through hardening and companion ex- 
changes to that observed there. This possibility needs to be 
investigated using high-precision N-hoAy computations of 
young star clusters. 

Given such reasonably-well quantified estimates of the 
distribution functions of orbital elements of the primordial 
binary population, the problem remains as to how these dis- 
tributions functions can be understood theoretically as a re- 
sult of the star-formation process. Fisher (2004) notes that 
distribution functions similar to the ones derived above, and 
in particular a wide mass-ratio distribution, a very wide pe- 
riod range, and a thermal eccentricity distribution, are ob- 
tained quite naturally from a turbulent molecular cloud (see 
also Burkert and Bodenheimer, 2000). 

3. THE PROPERTIES OF PRESTELLAR CORES. 

The gas that is just-about-ready to form stars arranges 
itself into denser structures often called prestellar 'cores' 
(e.g., Myers and Benson, 1983). Often, a 'typical' prestel- 
lar core is described as having a radius '--^ 0.1 pc, density 
> 10^ cm^'^, and velocity dispersion ~ 0.5 km s^^. In 
fact though, the idea that such cores are 'typical' primar- 
ily arises from the relative ease with which nearby, isolated 
dense cores, that will each form fewer than a handful of 
stars, can be observed and modelled. It is in fact likely that 
accounting for the diversity in core properties is crucial to 
improving the match between theory and observations of 
the conversion of gas to (binary) stars. 

3.1 What is a 'core'? 
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In general, observations and theory have concentrated on 

isolated and coherent prestellar cores such as those found 
in low-mass star forming regions such as Taurus (due to the 
relative ease of observing such cores). It is not yet clear if 
legitimate analogs to these cores exist within the dense con- 
centrations of gas that form the clusters (e.g., Goodman et 
ah, 1998). In particular, the ~ 0.1 pc size of isolated cores 
would result in them having multiple dynamical encounters 
in the dense environment that forms 'typical' clusters such 
as Orion (e.g., Lada andLada, 2003). 

That we have not observed any 0.1 pc core analogues 
in dense clusters is not surprising. Even in very nearby 
clusters hke NGC1333 in Perseus (at ~ 300 pc), 0.1 pc 
is 1 arcmin, typical of single-dish resolution for tracers like 
NH3 and N2H"'", which map out gas with density > 10^ 
(e.g., Benson and Myers, 1989; Evans, 1999). Thus, to find 
meaningful dense structures on scales significantly less than 
1 arcmin, interferometry is required. Interferometers have 
definitively revealed sub-structure in the gas within clus- 
ters, but this substructure does not offer a one-to-one gas 
clump-to-star match the way observations of isolated cores 
do. Instead, regions forming many stars are associated with 
more dense gas than those forming fewer. This lack of one- 
to-one correspondence suggests that long-lived blobs asso- 
ciated with the formation of individual cores within clusters 
do not exist. 

Given this, is it reasonable to extend observations and 
simulations of isolated cores to more typical clustered star 
forming 'cores'? The answer is possibly. Whilst large iso- 
lated cores cannot exist in clustered star forming regions, 
the regions in which stars are thought to form are far smaller 
than the size of a whole core. Observationally, the size of 
binary systems is <a few hundred AU in agreement with 
theoretical expectations of the scale of fragmentation (see 
Section 4.1). Systems of this scale would be expected to 
interact on timescales of > Myr in a typical cluster which 
is significantly longer than the star formation timescale is 
thought to be (see the chapters by Di Francesco et al. and 
Ward-Thompson et al). So, while the details of the first 
stages of collapse in isolated cores are probably not appli- 
cable to most star formation, the details of the final stages 
of fragmentation and star formation occurring on few hun- 
dred AU scales quite possibly occur in relative isolation. 
However, continued accretion onto cores may significantiy 
effect the evolution of the iimer proto-stellar system de- 
pending on the details and timescale of core and star for- 
mation in clusters. 

With this in mind, we will continue to review the prop- 
erties of isolated pre-stellar cores. 

3.2 Rotation. 

Clearly, for cores to fragment some angular momentum 
must be present otherwise the cores will collapse onto a 
single, central point. The simplest source of this angular 
momentum is due to bulk rotation of the core. It is a rela- 
tively straightforward procedure to estimate the component 



of solid-body rotation present in a dense core by fitting for 

the gradient in observed line-of-sight velocity over the face 
of a core (e.g., Goodman et al., 1993; Barranco & Good- 
man, 1998; Casein et al, 2002). The results of tiiis fitting 
(see Fig. 4) have been used as input values of 'initial angu- 
lar momentum' in many calculations. While the estimates 
of the component of soUd-body rotation made in this way 
are sound, and are thus fine to use as inputs, it is important 
to appreciate that cores do not really rotate as solid-bodies 
(Burkert and Bodenheimer, 2000). When the velocity mea- 
surements are put into the context of measurements of ve- 
locities on larger scales, both observations (Schnee et al, 
2005) and simulations {Burkert and Bodenheimer 2000), 
show that the "rotation" is often just an artifact created by 
larger-scale turbulent motions. 

Fig. 4 gives a summary of the measured specific an- 
gular momentum for core rotation for the 29 dense cores 
from Goodman et al. (1993) which show significant rota- 
tion. The majority of cores included in the Goodman et al. 
(1993) study are isolated, low-mass cores: one should keep 
in mind that the rotational properties of smaller fragments 
that may form inside those cores as true(r) precursors to 
protostars remain largely unmeasured. 

3.3 Non-thermal line widths. 

It has been known for many years {Larson, 1981; My- 
ers, 1983; Solomon et al., 1987; see also Elmegreen and 
Scalo, 2004a and references therein) that the line widths in- 
side even the most quiescent of dense cores are more than 
thermal. The coldest isolated dense cores have gas tempera- 
tures of order 10 K, and dust temperatures as low as 6 K (see 
the chapters by Di Francesco et al. and Ward-Thompson et 
al.). A gas temperature of 10 K implies an H2 la velocity 
dispersion of only 0.2 km s~^. Observed dispersions have a 
distribution from ~ 0.2 to 1 km s~^ with a peak at ~ 0.4 km 
s~^, never quite reaching down to the thermal vidue (e.g., 
see the catalogue of Jijina et at., 1999). 

The origin of the non-thermal Une width in dense cores is 
the subject of an extensive literature, but it is fair to say that 
a consensus exists that 'turbulence' is responsible (see e.g., 
Myers, 1983; Barranco & Goodman, 1998; Goodman et 
al., 199S; Elmegreen and Scalo, 2004-a). Significant levels 
of turbulence in cores are important for fragmentation as 
they may provide the angular momentum to form multiple 
systems (see Section 4.2.2) and, as mentioned above, could 
also be responsible for the observed rotation. 

3.4 Magnetic fields. 

The role of magnetic fields in supporting cores against 
collapse is a subject of much debate. Magnetic fields are 
not thought to be dynamically dominant in cores as was 
once thought (e.g., Shu et al., 1987; or indeed in the ISM as 
a whole Elmegreen and Scalo, 2004a). However, they may 
be very important as the role of magnetic fields in the frag- 
mentation of cores is poorly understood (see Section 4.2.6). 
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Isolated cores are found to be statistically triaxial with a 
tendency towards being prolate {Jones et al, 2001; Good- 
win et al., 2002). In contrast to this, magnetic support 
would tend to produce oblate cores. In addition these 
magnetically-supported oblate cores would tend to rotate 
around their short axis which is not observed (Goodman et 
al, 1993). This is supported by observations of the mag- 
netic field which show that cores are not magnetically criti- 
cal (Crutcher et al, 1999; see Bourke and Goodman, 2004 
and references therein). 

4. THE FRAGMENTATION OF PRESTELLAR 
CORES. 

4.1 The Physics of Collapsing Cores. 

During the early stages of the collapse of a prestellar 
core, the rate of compressional heating is low and the gas 
is able to cool radiatively, either by molecular line emis- 
sion, or, when p > 10^^^ g cm~^, by thermally coupling to 
the dust. The gas is therefore approximately isothermal (at 
~ 10 K) with an equation of state P oc p. 

Eventually the rate of compressional heating becomes so 
high (due to the acceleration of the collapse), and the rate of 
radiative cooUng so low (due to the increasing column den- 
sity and dust optical depth), that the gas switches to being 
approximately adiabatic, with P (x p^ (where 7 = 5/3 
initially for a monatomic gas, and then 7 = 7/5 above 
~ 300 K when H2 becomes rotationally excited). 

This behaviour has been studied in detail (Larson, 1969; 
Tohline, 1982; Masunaga et al, 1998; Masunaga andlnut- 
suka, 2000) for cores in the range 1 — lOA'/0 with an initial 
temperature of 10 K. These authors find that the switch be- 
tween isothermality and adiabaticity occurs at a critical den- 
sity of pcrit ~ 10"^^ g cm-3. (See Fig. 1 in Bate 1998 or 
Fig. 2 in Tohline, 2002). Thus, as contraction proceeds and 
the density, p, increases the Jeans mass, Mj oc p^^/^T^/^, 
decreases as long as the gas can retain the same tempera- 
ture, T, while it is optically thin. Once the opacity increases 
such that the gas core heats up, Mj increases. The most 
important result of this thermal behaviour is therefore that 
there is a minimum Jeans mass that is reached at ~ pcrit of 
order Mmin ~ 1O~^M0. This is often referred to as the 
opacity limit for fragmentation. 

There is an even lower minimum mass that occurs dur- 
ing a later isothermal phase at p ~ 10~^ g cm~^ when 
molecular hydrogen dissociates at a few thousand K. It is 
possible that a further fragmentation episode can occur at 
these densities which may account for some close binaries. 

Fragmentation in cores is expected to occur at around 
Pcrit as at lower or higher densities the Jeans mass increases 
(although how rapidly it rises above Pcr\% does depend sen- 
sitively on the 7 used in the adiabatic equation of state). 
Thus we expect multiple systems to be formed with a typi- 



cal length scale Rtorm of 

iiw < (i^)'^' - 125 (Meore/M0)l/3 AU (4) 

where Afporc is the mass of the core. Interestingly, this scale 
matches the observed peak in the T Tauri separation distri- 
bution (see also Sterzik et at., 2003). 

There is a minimum separation in this picture of ~ 
30 AU which is the separation of two fragments of Mmin 
at pcrit- It may - or may not - be significant that this is the 
average binary separation (DM91; see also Sterzik et al, 
2003). However, it would appear difficult to form binaries 
closer than ~ 20 — 30 AU without some hardening mecha- 
nism or a secondary fragmentation phase. 

It should be noted that the length scales of star formation 
of less than a few hundred AU are several orders of magni- 
tude smaller than the thousands of AU scales on which core 
properties have been observed. 

4.2 Fragmentation mechanisms. 

In this section we examine the main mechanisms that 
have been proposed to explain multiple formation. Given 
the complex and highly non-linear nature of the physics 
in most models, numerical simulations are the main route 
by which the mechanisms for fragmentation have been in- 
vestigated. Bulk rotation and turbulence are the two main 
mechanisms that have been considered to provide the angu- 
lar momentum required for fragmentation to occur and we 
review the theoretical work and simulations conducted on 
both of these mechanisms. In addition, we discuss the pos- 
sible role of magnetic fields, disc fragmentation and 'sec- 
ondary fragmentation' . 

4.2.7 Rotational fragmentation. The simplest situation 
in which fragmentation may well occur is in a spherical 
cloud with sohd-body rotation and an isothermal equation 
of state. Tohline (1981), using semi-analytic arguments 
concluded that all such clouds should fragment. A num- 
ber of simulations have shown that such clouds do fragment 
if athermAot :$0.12 - 0.15, where -£^therm/|f^| 
is the initial thermal virial ratio (where iJtherm is the ther- 
mal kinetic and O is the gravitational potential energy), and 
/3rot = £'rot/|f^| the initial rotational virial ratio (where 
i^rot is the rotational kinetic energy) (Miyama et al, 1984; 
Hachisu and Eriguchi, 1984, 1985; Miyama, 1992; see also 
Tsuribe and Inutsuka, 1999a, b; Tohline, 2002). 

Boss and Bodenheimer (1979) added an m = 2 az- 
imuthal density perturbation to a standard rotating cloud 
(effectively creating an elongated cloud more similar to 
those observed than purely spherical clouds; see Sec- 
tion 3.4). They found that with a perturbation of ampli- 
tude A = 0.5 the cloud fragments into a binary system. 
This simulation was repeated by Burkert and Bodenheimer 
(1993) who also found that when A = 0.1 a filament con- 
necting the two components of the binary fragments into 
several smaller fragments. However the connecting fila- 
ment should not fragment as predicted by Inutsuka and 
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Miyama (1992) and demonstrated by Truelove etal. (1997). 
Indeed the 'Boss and Bodenheimer test' has become a stan- 
dard test for the accuracy of codes (e.g., Truelove et al., 
1997 for adaptive mesh refinement - AMR - and Kitsonias 
and Whitworth, 2002 for smoothed particle hydrodynamics 
- SPH). However, it is a rather unsatisfactory test as, whilst 
the Truelove et al. simulations are generally considered to 
have converged, no analytic solution to the problem exists. 
An alternative test based on the original analysis of Jeans is 
presented by //wZ^foer a/. (2006). 

The simulation of rotating clouds can be made more 
physical by including an adiabatic (e.g., Tohline, 1981; 
Miyama, 1992) or barotropic (e.g., Bonnell, 1994; Bate and 
Burkert, 1997; Boss et al, 2000; Cha and Whitworth, 2003) 
equation of state (eos). Bate and Burkert (1997) showed 
that the Boss and Bodenheimer test does produce a line 
of fragments with a barotropic eos, but not if it remains 
isothermal. In addition, Bois ef a/. (2000) simulated a cloud 
with an TO = 2, A = 0.1 perturbation using a barotropic eos 
and also with radiation transport; the second case producing 
a binary whilst the first did not: despite the similarity of the 
pressure-temperature relations. Both of these results sug- 
gest that fragmentation is highly sensitive to thermal inertia 
and radiation transport effects. 

Other authors have modified the initial conditions to in- 
clude effects such as different density profiles (e.g., Myhill 
and Kaula, 1992; Burkert et at., 1997; Boss, 1996; Boss 
and Myhill, 1995; Burkert and Bodenheimer, 1996; Boss et 
al, 2000; Boss, 1993), differential rotation (which tends to 
promote fragmentation: Myhill and Kaula, 1992; Boss and 
Myhill, 1995; Cha and Whitworth, 2003) and non-spherical 
shapes (e.g., Bastien, 1983; Bonnell and Bastien, I99I; 
Bonnell et al., 1991; Nelson and Papaloizou, 1993; Boss, 
1993; Sigalotti and Klapp, 1997). The effect increasing ex- 
ternal pressure on the collapse of rotating cores have been 
investigated by Hennebelle et al (2003, 2004, 2006). 

4.2.2 Turbulent fragmentation. Recently a picture of 
star formation as a rapid and highly dynamic process has 
appeared (e.g., Elmegreen, 2000; Vdquez-Semadeni et al, 
2000; Larson, 2003; Elmegreen and Scalo, 2004b) as op- 
posed to a quasi-static process (e.g., Shu et al, 1987). In 
particular, the idea of cores evolving slowly via ambipo- 
lar diffusion (e.g., Basu and Mouschovias, 1994, 1995a, 
b; Ciolek and Mouschovias, 1993, 1994, 1995; Ciolek 
and Basu, 2000) has been replaced by one in which cores 
form in converging flows in a highly turbulent molecular 
cloud. This is rather good news for fragmentation, as the 
main effects of a quasi- static evolution are to delay frag- 
mentation and reduce the angular momentum and turbu- 
lence in a core and organise material so that its coUapse 
is well focused onto a central point (see also Section 3.2). 
Simulations of core formation in a turbulent medium sug- 
gest that cores form with significant amounts of turbulence. 
Turbulent, rapidly formed cores also reproduce many of 
the observed properties of cores (e.g., Burkert and Boden- 
heimer, 2000; Ballesteros-Paredes et al , 2003 ; Jappsen and 
Klessen, 2004) and have a mass spectrum not dissimilar to 



the observed core mass spectrum (e.g., Padoan and Nord- 
lund, 2002, 2004; Klessen et al, 2005). 

Simulations of the effects of turbulence in cores fo- 
cus on two different regimes: high- velocity (Bate et al, 
2002, 2003; Bate and Bonnell, 2005; Delgado Donate et 
al, 2004a, b), and low- velocity (Goodwin et al., 2004a, 
b). Also see Fisher (2004) for a semi-analytic approach 
to multiple formation with turbulence. The level of turbu- 
lence is usually quantified as a turbulent virial ratio aturb = 
Etuib/\^\< where .Bturb is the kinetic energy in turbulent 
motions and is the gravitational potential energy (note - 
not any rotational property). Highly turbulent simulations 
focus on aturb = 1 in 50Mq (Bate et al., 2002, 2003; 
Bate and Bonnell, 2005) and 5Mq (Delgado Donate et al, 
2004a, b) cores. Simulations of slightly turbulent cores 
range between aturb = — 0.25 in 5.4M0 cores (Good- 
win et al, 2004a, b). In all of these simulations turbulent 
motions are modelled using a Gaussian divergence-free ran- 
dom velocity field P{k) oc fc^" where n is usually taken to 
be 4 to match observations of cores for which n = 3—4 pro- 
vides a good fit to the Larson relations (Burkert and Boden- 
heimer, 2000). It should be noted that the random chaotic 
effects introduced by variations in the initial turbulent ve- 
locity field can be very important. Therefore a statistical 
approach is desirable utilising large ensembles of simula- 
tions (e.g., Larson, 2002). 

These simulations of turbulence are different to those of 
turbulence in molecular clouds which concentrate on the 
formation of dense cores and massive stars (e.g., Klessen 
and Burkert, 2000). This is largely due to computational 
limitations which do not allow the resolution of the opacity 
limit for fragmentation in the large-scale context of giant 
molecular clouds. A mass resolution of ~ 1O~^M0 is re- 
quired to resolve the opacity Umit for fragmentation in SPH. 
In AMR the problem is even worse as the Jeans length con- 
tinues to fall (albeit more slowly) after the minimum Jeans 
mass is reached and codes must resolve few AU scales to 
capture the lowest-mass fragments. 

Even very low levels of turbulence (aturb ~ 0.025) are 
enough to allow cores to fragment: that is, for most cores 
in an ensemble to form more than one star (Goodwin et al , 
2004b). As the level of turbulence is increased, the average 
number of stars that form in a core increases (Goodwin et 
al, 2004b). It has been suggested that approximately one 
star forms per initial Jeans mass: ~ IMq for these ini- 
tial temperatures and densities (cf. Bate et al, 2003; Del- 
gado Donate et al, 2004a). This seems to hold in highly- 
turbulent cores, however the number of stars forming falls 
with decreasing turbulence (Goodwin et al, 2004b) and so 
this at best probably only represents a (statistical) asymp- 
totic behaviour. 

In highly-turbulent cores, the supersonic turbulent ve- 
locity field creates a number of condensations in shocked, 
converging regions which become Jeans unstable and col- 
lapse (see Bate et al, 2003; Delgado Donate et al, 2004a, 
b). However, it is unclear if this mode of fragmentation 
is realistic in small (certainly < a few Mq) cores as the 
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observed levels of non-themal motions rule-out significant 
highly supersonic turbulence in these cores. 

In Fig. 5 we show the formation of a fragment in a 
mildly-turbulent 5AMq core with Ofturb = 0.05 based on 
observations of the isolated core L1544 (from Goodwin et 
al, 2004a). Fragmentation occurs in a 'disc-like' mode in 
circumsteUar accretion regions (we avoid the use of 'disc' 
to describe these regions as they are not rotationally sup- 
ported structures) (CARs) which form around the first star. 
CARs are highly unstable structures as there is non-uniform 
(in space, time and angular momentum) inflow onto them. 
Complex spiral instabilities form in the CAR due to this 
inhomogeneous infall of material. We note that these in- 
stabilities are seen in both SPH and AMR simulations of 
the same situation (Gawryszczak et al, 2006; see also Sec- 
tion 4.4). Fragmentation occurs if the density in spiral 
waves becomes high enough that the Jeans length falls to 
the typical width of a spiral wave and the collapse time falls 
to a low enough fraction of the local rotation period that 
it may escape shredding by differential rotation. In these 
simulations, fragmentation occurs for some (but not all) re- 
gions that exceed ^ 10^^^ g cm^'^ in density (equating to 
a Jeans length of 20 AU) beyond 50 - 100 AU from 
the central star. We note that the highly unstable nature of 
CARs makes usual applications of instabiUty criteria such 
as the Toomre Q-parameter impossible. 

Such a mode of fragmentation is highly sensitive to the 
equation of state that has been adopted. It has been found 
that the number of fragments that form increases if 7 is 
changed from 5/3 to 7/5 {Goodwin et al, in prep.). This is 
due to the sensitivity of the Jeans length with density and so 
to the ease with which fragmentation can occur in CARs. 

The process of fragmentation in CARs is highly chaotic, 
relying as it does on a certain degree of 'luck' in being 
able to reach a high-enough density and avoiding shredding 
whilst collapsing. Thus it is no surprise that anywhere be- 
tween 1 and 12 stars form in each core depending entirely 
on the details of the initial turbulent velocity field {Goodwin 

a/., 2004a, b). 

In sunmiary, it is found that turbulent cores generally 
fragment into several stars: approximately one per initial 
Jeans mass (~ IMq) in the core. The number of stars 
that form increases with increasing turbulence and is also 
highly sensitive to the details of the turbulent velocity field. 
However, only relatively high-mass cores (> 5Mq) have 
been investigated in turbulent simulations so far. The effect 
of turbulence in lower-mass cores must be investigated, as 
lower-mass cores appear to dominate the core mass func- 
tion {Motte et al, 1998; Testi and Sargent, 1998; Matte et 
al, 2001). 

4.2.3 Disc fragmentation. Disc fragmentation is a mech- 
anism by which low-mass stars and BDs may be formed. 
In the dense envirormients of clusters close encounters be- 
tween stars can disturb the circumsteUar discs promoting 
instabilities which can lead to the fragmentation of oth- 
erwise stable discs. (Note that this is rather different to 
the turbulent disc-like scenario described above as these 



proto-planetary discs are much less massive than CARs and 
are also stable, rotationally supported discs as opposed to 
CARs). 

In a series of papers. Boffin et al. (1998) and Watkins 
(1998a,b) found that most star-disc interactions will lead to 
gravitational instabilities which form new low-mass com- 
panions. These simulations generally considered massive 
discs where M^tar = ^^disc = O.5M0. Bate et al. (2003) 
find that star-disc encounters play an important role in form- 
ing binaries and also truncating discs. Star-disc encounters 
are also thought to play an important role in redistributing 
angular momentum in proto-planetary discs even if they do 
not cause further fragmentation {Larson, 2002; Pfalzner, 
2004; Pfalzner et al, 2005). 

Star-disc encounters probably play a role in star forma- 
tion, and may lead to to the formation of BD (or even plan- 
etary) mass companions {Thies et al, 2005). However, they 
are probably not a significant contributor to the primordial 
stellar binary population. This is due to the requirement 
that the encounters occur early in the star formation process 
- during the Class phase when the disc mass is still very 
large compared to the stellar mass - a phase which lasts for 
only ^ 10'' yrs, leaving only a small time for encounters 
to occur However, the role of disc fragmentation in planet 
formation may well be important. 

4.2.4 The role of magnetic fields. The treatments of col- 
lapse and fragmentation discussed above do not include 
magnetic fields. The new picture of rapid, turbulence- 
driven star formation combined with the lack of observa- 
tional evidence for magnetically critical cores suggests that 
magnetic fields are not dynamically dominant. In addition, 
fragmentation is expected to occur at densities ^ 10~^^ g 
cm~^, densities at which the magnetic field is expected to 
be decoupled from the gas due to the extremely low frac- 
tional ionisation (see Tohline, 2002). However, possibly 
one of the main reasons for neglecting magnetic fields is the 
difficulty in including them in SPH simulations, (although 
this is improving, see esp. Hosking and Wliitworth, 2004a, 
b; Price and Monaghan, 2004a, b). Magnetic fields are 
clearly present in (many) cores, even if they are not dynam- 
ically dominant, and their effects may be very important. 

Grid-based simulations which include magnetic fields in 
rotating clouds show that fragmentation can occur in these 
clouds, although magnetic fields appear to have a tendency 
to suppress fragmentation (e.g., Hosking and Whitworth, 
2004b; Machida et al, 2005b), although Boss (2002;2004) 
claims the opposite. Sigalotti and Klapp (2000) find bi- 
nary and higher-order multiple formation in slowly rotating 
^ IMq clouds which includes a model for ambipolar dif- 
fusion. 

Possibly the most extensive investigation of the effects 
of magnetic fields on fragmentation has been made by 
Machida et al (2005a,b). They find that fragmentation oc- 
curs in ~ 50% of their rotating, magnetised clouds when ei- 
ther the rotation is relatively high or magnetic field sttength 
relatively low. In particular fragmentation always occurs 
in magnetised clouds if /3rot > 05, but it almost never 
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Fig. 5. — The evolution of a CAR in a turbulent molecular core from Goodwin et al. (2004a). The boxes are 1000 AU on a side and the 
timescale runs from 66.6 to 67.8 kyr after the start of the simulation in steps of 400 yrs. The grey-scale bar gives the column density in 
g cm^^. The spiral features can become self-gravitating if their density exceeds ~ 10"^^ g cm~"^ as their Jeans length falls to ~ 20 AU 
which can allow collapse without being shredded. This does not always occur, in the first two panels a dense knot can be seen being 
shredded and accreted onto the central object. Very similar evolution is also seen in AMR simulations {Gawryszczak et al, 2006). 



occurs below this limit (see Fig. 10 from Machida et al, 
2005b). Indeed, Burkert and Balsam (2001) conclude that 
once magnetic fields are strong enough to affect the dynam- 
ical evolution they will also efficiently suppress fragmenta- 
tion which means that magnetic fields cannot be important 
as we know that fragmentation must occur. 

4.3 'Secondary' fragmentation. 

As briefly mentioned in Section 4.1, there is a second 
isothermal phase in the evolution of gas towards steUar den- 
sities. This occurs at a temperature of 2000 K and a 
density of ^ 10~^ g cm~^ when molecular hydrogen dis- 
sociates into atomic hydrogen. This phase occurs in the 
hydrostatic protostar when its radius is ^ 1 AU and - if 
fragmentation can occur at this stage - it may explain very 
close binaries. 

Bofli Boss (1989) and Bonnell and Bate (1994) simu- 
lated the collapse of a rotating hydrostatic first object to 
high densities. They found that fragmentation can occur in 
axisymmetric instabiUties or a ring formed by a centrifugal 
bounce. However, Bate (1998) found that spiral instabili- 
ties remove angular momentum and suppress further frag- 
mentation. Recent 2D simulations by Saigo and Tomisaka 
(2006) suggest that the angular momentum of the first core 
is a crucial factor in determining if fragmentation will occur 
during the second collapse. 

Thus it is unclear if a secondary fragmentation phase oc- 
curs. However, we suggest that such a phase could well be 
responsible for the apparently high incidence of very close 
BD-BD binary systems (Pinfield et al, 2003; Maxted and 
Jeffries, 2005) as the evolution of BD-mass hydrostatic ob- 
jects occurs on a longer timescale than in stellar-mass ob- 
jects. This possibility is being investigated by SG without 
any firm conclusions as yet. 

4.4 Simulations vs. Observations 

A sunamary of the simulations to date suggests that col- 
lapsing cores are easily able to fragment. However, no de- 



tailed model is currently able to correctly predict all of the 
observed binary properties. 

A successful model of star formation must produce mul- 
tiple systems which generally have only 2 or 3 stars with 
a wide range of separations from << 1 AU to a peak at 
~ 100 — 200 AU. At all separations, most stars must usu- 
ally have quite different masses, but avoiding BDs within at 
least 5 AU of the primary (the BD desert). 

Possibly the most significant problem at the moment, is 
that simulations seem to form too many single stars (see 
Bouvier et al, 2001; Duchene et al, 2004; Goodwin and 
Kroupa, 2005). As described in Section 2.3.1 systems with 
iV > 3 are generally unstable and decay by ejecting their 
lowest mass member and hardening the remaining multi- 
ple. The ejection of members of small- A/^ multiple systems 
dilutes the multiplicity of stars, as ejected stars tend to be 
single. Thus, many ejections will result in a far lower mul- 
tiplicity fraction than is observed in young star forming re- 
gions. Goodwin and Kroupa (2005) suggest that the ob- 
served multipUcity frequencies can be explained if roughly 
half of cores form 2 stars, and half form 3 stars. However, 
these numbers are far lower than are usually found in core 
fragmentation simulations. 

The inclusion of magnetic fields produces the opposite 
problem that too few binaries are produced. Machida et 
al (2005b) find the fragmentation does not occur in rotat- 
ing, magnetised clouds when /3rot ~ 0.04 - a higher level 
of rotation than is observed in many cores. This problem 
becomes especially accuse when we consider that much of 
the observed rotation in cores could well be due to turbulent 
motions rather than a bulk rotation. 

A related problem is posed by the existence of a signif- 
icant number (~ 20% of field G-dwarfs; DM91) of close, 
unequal mass binary systems. It appears difficult to form 
stars much closer together than ~ 30 AU. In order to ob- 
tain hard binary systems, a further hardening mechanism 
is required. In many simulations of turbulent star forma- 
tion, this hardening mechanism is provided by the ejection 
of low-mass components (e.g.. Bate et al., 2003; Delgado 
Donate et al., 2004a, b; Goodwin et al., 2004a, b; Um- 
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Fig. 6. — Comparison of an AMR (left) and a SPH (right) simulation of a collapsing, turbulent 5AMq core showing the first 'disc' 
fragmentation episode from Gawryszczak et al. (2006). The bound fragments can be seen at (0,125) AU in the AMR simulation (left) 
and (-200,0) AU in the SPH simulation (right). The scale and orientation of both views are identical. 



breit et al, 2005; Hubberand Whitworth, 2005; see Ochi et 
al, 2005 for a caveat). However, ejections appear to occur 
rapidly, during the main accretion phase, producing many 
close, equal-mass binary systems (see above) in contradic- 
tion to observations of a relatively flat mass ratio distribu- 
tion {Mazeh et al, 1992). White and Ghez (2001) do find 
many roughly equal-mass PMS binaries < 100 au in Tau- 
rus, however there is no trend to more equal mass ratios at 
very low separations. Fisher et al (2005) do find a bias to- 
wards equal-mass binaries among local field spectroscopic 
binaries (with separations < 1 AU). 

In general, simulations that produce a small number of 
stars consistent with limits on ejection do not produce hard 
binaries. These binaries are difficult to form in significant 
numbers through later dynamical interactions in a clustered 
environment, which tend to disrupt wide binaries but not 
harden them. But simulations that produce many stars tend 
to form too many close, equal-mass binaries and very high- 
order multiple systems (many quadruples and quintuples) 
and also dilute the multiplicity fraction too much through 
ejections. 

4.5 Numerical issues 

As already discussed in this volume by Klein et al , there 
is some debate about the ability of simulations to correctly 
resolve fragmentation. Simulations of star formation are 
usually conducted using SPH as opposed to AMR schemes 
due to the Lagrangian nature of SPH (see Gawryszczak et 
a/., 2006 for details). 

No numerical scheme is perfect, and both SPH and AMR 
have their advantages and disadvantages. It is worth not- 
ing that a recent study by Gawryszczak et al (2006) has 
shown that AMR and SPH converge when simulating the 
collapse of a slightly turbulent core. Fig. 6 shows a snapshot 
of both the SPH and AMR simulations at roughly 76 kyr 
from the start of the simulation showing that in both numer- 



ical schemes a highly unstable CAR forms that fragments 
in both simulations. Gawryszczak et al found that AMR 
is significantly more computationally intensive than SPH 
for an identical simulation. This result is not surprising as 
when simulating gravitational collapse the Lagrangian na- 
ture of SPH should prove highly efficient. The agreement of 
two very different methods when applied to the same physi- 
cal situation increases our confidence that the results are not 
dominated but numerical effects. 

Both SPH and AMR suffer from problems with artificial 
angular momentum transport. In AMR, rotation in a poorly 
resolved Cartesian grid is likely to transport angular mo- 
mentum outwards. In SPH, the use of artificial viscosity to 
reduce particle inter-penetiation in shocks produces an in- 
ward transport of angular momentum with rotation. These 
problems are probably responsible for the different CAR 
(disc) sizes between SPH and AMR seen by Gawryszczak 
et al (2006; see also Fig. 6). 

It should be noted that Hubber et al (2006) have shown 
that SPH supresses artificial fragmentation rather than pro- 
moting it which suggests that the current generation of SPH 
simulations could underestimate the number of fragments 
that form. 

We would conclude that the computational situation is 
not perfect and many problems with both SPH and AMR 
remain. However, we think that the conflict between sim- 
ulation and observation is more probably due to missing 
and/or incorrect physics, rather than any fundamental nu- 
merical difficulties. 

The computational situation is made more problematic 
by the need to perform ensembles of simulations to get the 
statistical properties of multiple systems. Even in cases 
where the result for any single set of parameters might be 
expected to converge there is a large parameter space to 
cover, and small differences in initial conditions may make 
a significant difference to the result. More reahstically, sit- 
uations where fragmentation occurs due to non-Unear insta- 
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bilities, or with turbulence (when the details of the velocity 
field vary) require large ensembles even for the same re- 
gion of parameter space. No matter how detailed or correct 
any single simulation is, it can only be a snapshot of the 
outcome of a particular initial configuration. This vastly 
increases the computational effort required to model frag- 
mentation and star formation. 

4.6 Missing physics 

One of the greatest problems facing the simulation of 
core fragmentation is to correctly model the thermal physics 
of cores. It is not possible in the foreseeable future that we 
will be able to conduct hydrodynamical simulations that in- 
clude a proper treatment of radiation transport as the com- 
putational expense is just too great. However, as shown by 
Boss et al. (2000) thermal effects can be very important. 
The use of the barotropic equation of state is at best a first 
approximation, and it is clear that varying the adiabatic ex- 
ponent can have significant effects on fragmentation {Good- 
win et al, in prep.). In particular, the barotopic equation 
of state is based on simulations that have used (necessar- 
ily) simplistic, spherically symmetric assumptions. It is not 
clear to what extent these can be applied to highly inhomo- 
geneous cores including discs and local density peaks. Im- 
provements are being made, including using a flux-limited 
diffusion approximation (Whitehouse and Bate, 2004) and 
an approximation based on the local potential as a guide to 
optical depth (Stamatellos et al., in prep.). However a prob- 
lem remains with any approximation in that, while it should 
obviously match the fully detailed radiative transfer simu- 
lations of simple situations, it is not clear if it is correct in 
more complex situations. 

Most simulations (especially SPH simulations) do not 
include the effects of magnetic fields. Yet those simula- 
tions which do include magnetic fields seem to suggest that 
fragmentation is suppressed. This could well be a very im- 
portant conclusion given that non-magnetic models seem 
to over-produce stars in cores. However, the very efficient 
suppression of fragmentation by magnetic fields may rule- 
out the importance of magnetic fields in the fragmentation 
process as we know that cores must fragment. 

Given that fragmentation appears to occur in disc-like 
structures, the proper treatment of these is vital. Both AMR 
and SPH have problems with the artificial transport of an- 
gular momentum (see Gawryszczak et al, 2006) which will 
effect their ability to correctly model discs. Discs are also 
a situation in which magnetic fields may play an important 
role. 

Finally, very few simulations attempt to model the ef- 
fects of feedback from stars as jets or through their radiation 
field (for some first attempts to deal with these problems see 
Stamatellos et al, 2005; Dale et al, 2005). In particular, 
Stamatellos et al. (2005) find that the inclusion of jets may 
inhibit fragmentation by decreasing the inflow rate onto the 
disc and forcing that inflow to occur away from the poles. 

Many of the physical situations which may result in frag- 



mentation are rather complex and often chaotic (turbulence 

being the most obvious example). Such situations will not 
produce any single, unique answer. Indeed, given the va- 
riety of multiple systems such a situation would not be 
expected. However, this does require that a statistical ap- 
proach be taken when performing simulations. This vastly 
increases the computational effort required, as any 'single' 
region of an already huge parameter space will require an 
ensemble of simulations to investigate it. 

5. CONCLUSIONS 

Almost all young stars are found in multiple systems 
with a very wide separation distribution and a fairly flat 
mass ratio distribution. Thus prestellar cores must fragment 
into multiple stars and/or BDs with these properties. 

The dynamical decay of small- A?^ systems would rapidly 
produce a large single-star pre-main sequence population if 
large numbers of unstable systems form with A'" > 2 or 3. 
This decay would also result in large numbers of very close 
binary systems. Neither of these are observed, leading to 
the conclusion that cores must usually form only 2 to 4 stars 
in hierarchical systems (for N > 2). 

Simulations show that most cores which contain some 
angular momentum - either in bulk rotation, or in turbu- 
lence - are able to fragment into multiple objects. However, 
these simulations have been unsuccessful in matching their 
results to the observed young multiple population. In par- 
ticular, the distributions of separations and mass ratios from 
simulations tend not to fit well. 

The future is somewhat rosier, however. The inclusion of 
more detailed physics and more realistic initial conditions 
may well yield better fits to observations. 
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